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MODELING  THE  DYNAMIC  LOAD/UNLOAD  BEHAVIOR 
OF  CERAMICS  UNDER  IMPACT  LOADING 


J.R.  Furlong,  J.F.  Davis,  and  M.L.  Alme 
R  &  D  Associates 

ABSTRACT: 

The  RDA  ceramic  failure  model  was  used  to  model  the  behavior  of  SiC  and  B4C 
during  plate  impact  loading.  The  experiments  modeled  were  performed  at  Sandia  National 
Laboratory  (Kipp  and  Grady,  1989).  The  simulations  were  performed  with  a  version  of  the 
Sandia  WONDY-V  hydrocode  (Kipp  and  Lawrence,  1981)  that  was  modified  to  run  on  a  PC 
and  to  include  the  RDA  ceramic  failure  model. 

The  model  adequately  characterizes  the  two  ceramic  materials  although  B4C  and  SiC 
exhibit  distinctly  different  degrees  of  work  hardening,  strain  rate  hardening,  and  damage 
sensitivity.  Both  materials  experienced  damage  (micro-cracking)  on  load  which  was  predicted 
by  the  model.  The  experimental  unloading  data  exhibits  a  reduced  magnitude  elastic  unload 
wave  which  suggests  that  further  damage  has  accumulated  after  the  ceramic  reaches  peak 
stress.  A  strength  reducing  mechanism  such  as  fatigue  or  melting  at  grain  boundaries  is 
suspected  as  the  cause  of  this  unpredicted  strength  loss. 

No  spall  was  evident  in  the  measured  signals.  The  model  predicts  spall  will  not  occur 
even  though  the  predicted  axial  stress  levels  were  significantly  above  the  estimated  spall 
strengths  of  the  two  ceramics.  Failure  was  prevented  by  positive  (compressive)  mean  stress 
levels  that  were  maintained  throughout  the  ceramic  specimen  during  the  experimental 
measurement.  A  suggested  conclusion  is  that  for  tensile  failure  to  occur  in  these  ceramics, 
the  mean  stress  must  be  negative  (tensile). 
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L  INTRODUCTION: 


High-strength  ceramics  have  been  shown  to  be  effective  armor  materials  and  therefore, 
the  characterization  of  their  dynamic  behavior  is  of  great  interest.  The  RDA  ceramic  failure 
model  (Furlong,  Alme  and  Rajendran,  1988)  has  been  shown  to  be  a  valuable  tool  in 
simulating  ceramics  under  ballistic  impact  loading  and  developing  design  parameters  for 
armor  system  design.  The  modeling  of  the  dynamic  load/unload  behavior  of  ceramics  under 
plate  impact  loading  were  performed  with  a  version  of  the  Sandia  WONDY-V  hydrocode 
(Kipp  and  Lawrence,  1981)  that  was  modified  to  run  on  a  PC  and  to  include  the  RDA 
ceramic  failure  model.  The  major  focus  of  the  modeling  program  has  been  to  model  SiC  and 
B4C  plate  impact  experiments  performed  at  Sandia  National  Laboratory  (Kipp  and  Grady,  1989). 

A  systematic  approach  to  modeling  these  materials  was  adopted  which  entails 
assembly  of  available  Hugoniot  data,  evaluation  of  the  data  to  identify  a)  the  underlying 
hydrostat,  b)  any  unusual  strength  behavior,  c)  the  effect  of  porosity  on  the  hydrostat,  and 
d)  any  evidence  of  phase  transformations.  Next,  with  an  assumed  hydrostat,  the  material  was 
modeled  as  an  elastic  perfectly-plastic  material,  then,  noting  differences  between  the  elastic 
perfectly-plastic  solution  and  the  observed  particle  velocity  history,  make  adjustments  to  the 
RDA  model  constants  for  strength  and  damage  until  a  suitable  match  was  made  (or  could  not 
be  made).  Finally,  incorporate  improvements  into  the  model  to  reflect  new  insights  of 
ceramic  material  properties. 

Details  of  this  analysis  procedure  are  described  in  separate  sections  of  this  report. 

Worthy  of  note  was  the  unexpected  result  that  none  of  the  ceramics  modeled  spalled  during 
the  period  of  measurement.  This  result,  predicted  by  the  model,  was  independent  of  the 
impact  velocities  studied,  and  occurred  at  axial  tensile  stress  levels  as  high  as  25  kbar  (well 
beyond  the  3-5  kbar  range  normally  associated  with  spall  strength  of  low  porosity,  high 
strength  ceramics).  Common  in  all  the  simulations  was  the  maintenance  of  a  positive 
(compressive)  mean  stress  (pressure)  throughout  the  specimen  that  was  primarily  the  result  of 
the  close  impedance  match  between  the  ceramic  tiles  and  the  lithium  fluoride  window  of  the 
VISAR  instrumentation.  Further  analysis  of  the  model  results  indicates  that  for  tensile 
failure  in  these  ceramics  to  occur,  the  mean  stress  must  be  negative  (tensile). 
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2.  REVIEW  OF  RDA  CERAMIC  FAILURE  MODEL: 


The  RDA  ceramics  failure  model  is  a  combined  tensile-compressive  failure  model. 

The  tensile  and  compressive  phases  each  have  their  own  set  of  constitutive  laws  and  damage 
evolution  equations.  The  criterion  for  branching  into  either  model  is  the  sign  of  the  mean 
stress  (pressure).  If  at  the  beginning  of  each  integration  cycle  the  mean  stress  is  positive  then 
the  logic  for  compression  is  executed.  If  the  mean  stress  in  negative  then  the  logic  for 
tension  is  executed. 

The  underlying  concept  of  the  tensile  model  is  that  of  a  degraded  modulus  theory.  To 
implement  the  theory  the  usual  linear  elastic  stress  -strain  equations  are  modified  as 


P  =  3K  (1  -  D)  6 

(1) 

Sy  =  2G  (1  -  D)  ey 

(2) 

where  P  is  the  pressure  (mean  stress),  Sy  and  ey  are  the  deviatoric  components  of  stress  and 
strain,  respectively,  <)>  is  the  volumetric  strain,  K  and  G  are  the  undamaged  bulk  and  shear 
modulus  of  the  material,  and  D  is  a  state  variable  describing  the  internal  damage  in  the 
material.  These  equations  govern  the  entire  elastic  and  plastic  stress  response.  Concepts  such 
as  equations-of-state  and  yield  surfaces  do  not  apply.  Hence,  the  accuracy  of  the  computed 
stress  levels  is  based  on  the  accuracy  of  the  damage  evolution  equation. 

The  tensile  damage  model  is  microphysically  based  and  adopts  a  self-consistent 
approach  to  describe  the  elastic  modulus  of  a  cracked  body  by  a  crack  density  (void  volume) 
that  is  a  function  of  a  degraded  Poisson’s  ratio  (Budiansky  and  O’Connell,  1976).  When  the 
crack  density  reaches  a  level  of  9/16,  or  slightly  more  than  one-half  the  total  volume  of  an 
elemental  section  of  material,  the  modulus  is  totally  degraded  and  the  element  ceases  to  have 
any  stiffness. 

The  RDA  compressive  damage  model  is  a  plasticity  based  phenomenological  model. 
Johnson  and  Cook’s  (1983)  expression  for  the  yield  surface  of  a  ductile  material  has  been 
extended  to  include  a  damage  softening  term,  (1-D).  The  resulting  expression  for  the 
damaged  softened  yield  surface  becomes; 

Y  =  ((Cj  +  C2  £pn)  (1+  C3  In  e*)  (1  -  T*m)  +  C4  P)  (1  -  D)  (3) 
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where  Cj  is  the  static  yield  strength,  C2  and  n  are  the  strain  hardening  parameters,  C3  is  the 
viscoplastic  (strain  rate)  parameter,  T*  and  m  are  the  thermal  softening  terms,  C4  is  the 
confining  pressure  parameter,  and  (1-D)  is  the  damage  induced  softening  term.  Therefore, 
this  model  can  be  used  to  calculate  the  strength  from  the  combined  effects  of  work  hardening, 
strain  rate,  thermal  softening,  and  confining  pressure.  Damage  evolution  is  assumed  to  be  a 
function  of  the  plastic  work  rate  in  compression  only.  A  single  damage  sensitivity  parameter 
is  used  to  characterize  any  given  materials’  dependence  of  damage  on  plastic  work  rate. 

Under  compressive  loading,  damage  affects  only  the  strength  of  the  material,  not  its  stiffness 
(no  change  in  the  longitudinal  modulus).  However,  the  combined  model  maintains  a  single 
array  for  damage,  so  that  damage  evolved  in  tension  will  affect  the  strength  of  the  material  in 
compression,  and  damaged  evolved  in  compression  will  affect  the  stiffness  of  the  material  in 
tension. 

A  more  complete  description  of  the  model  and  model  theory  is  included  as 
Appendix  1. 
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3.  MODELING  METHODOLOGY: 


A  systematic  approach  was  used  to  characterize  these  materials.  Thw  Initial  step  was 
to  assemble  all  available  Hugoniot  data,  and  to  evaluate  Hugoniot  data  for  identification  of 
a)  the  underlying  hydrostat,  b)  any  unusual  strength  behavior,  c)  the  effect  of  porosity  on  the 
hydrostat,  and  d)  any  evidence  of  phase  transformations.  The  next  step  was  to  model  the 
material  as  elastic  perfectly-plastic  with  an  assumed  hydrostat.  By  noting  the  differences 
between  the  elastic  perfectly-plastic  solution  and  the  observed  particle  velocity  history,  the 
adjustments  were  made  to  the  RDA  model  constants  for  strength  and  damage  until  a  suitable 
match  was  made  (or  could  not  be  made).  Finally,  improvements  into  the  model  to  reflect 
new  insights  of  ceramic  material  properties  were  incorporated  to  the  model. 

If  a  material  behaves  in  a  classic  elastic-perfectly  plastic  manner,  a  representative 
load-release  path  at  a  position  in  the  target  material  in  uniaxial  strain  would  have  the 
appearance  of  the  path  plotted  in  Figure  1.  Note  that  a  reference  hydrostat  is  also  included  in 
the  figure.  The  load  path  is  elastic  to  the  Hugoniot  elastic  limit  (HEL)  with  a  slope  equal  to 
the  longitudinal  modulus,  pCL2,  where  a  break  in  the  curve  occurs  at  the  HEL  and  the 
material  behaves  as  a  plastic,  and  loading  continues  as  defined  by  the  underlying  hydrostat, 
maintaining  a  separation  from  the  hydrostat  of  2/3  the  yield  strength.  When  unloading 
commences,  the  material  releases  elastically,  crosses  the  hydrostat,  and  continues  to  unload  to 
a  value  of  2/3  the  yield  strength  below  the  hydrostat  with  further  unloading  on  a  curve  as 
defined  by  the  hydrostat.  The  paths  are  characterized  by  sharp  transitions  from  elastic  to 
plastic  states.  As  discussed  in  the  previous  section,  the  RDA  ceramic  model  assumes  damage 
evolution  is  a  function  of  the  plastic  work  rate  in  compression.  Thus,  above  the  HEL  where 
the  material  behaves  as  a  plastic,  damage  to  the  ceramic  will  collapse  the  Hugoniot  to  the 
underlying  hydrostat.  This  effect  is  shown  schematically  in  Figure  1  with  the  dashed  line. 

Above  the  HEL,  the  relative  contributions  of  the  strength  and  pressure  is  not  known. 
The  assumption  must  be  made  of  the  form  for  the  hydrostat,  then  using  an  appropriate  yield 
function,  such  as  the  Johnson-Cook  model,  make  adjustments  with  the  strength  model  to 
predict  the  shape  and  velocity  of  the  plastic  wave.  For  many  materials,  the  shock  velocity  us 
may  be  expressed  as  a  linear  function  of  particle  velocity  u  as  given  by: 
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(4) 


us  =  C o  +  sup 

then 

P-PQ^I 

(1-sri)2 


(5) 


Therefore  the  hydrostat  can  be  based  on  a  simple  shock  velocity/particle  velocity  relation  with 
one  free  parameter,  s,  that  is  used  to  match  the  Hugoniot  data  at  high  compression  (where  the 
strength  of  the  material  is  assumed  negligible).  This  assumption  for  the  form  of  the  hydrostat 
assumes  a  linear  function  of  shock  velocity  versus  particle  velocity  over  the  stress-strain  space 
where  the  stress  is  the  sum  of  pressure  and  strength. 

Hugoniot  data  for  ceramic  materials  with  porosity  indicates  that  at  each  initial  value  of 
porosity,  the  ceramic  will  behave  as  a  different  material.  This  effect  is  shown  schematically 
in  Figure  2.  Depending  on  the  initial  conditions  of  the  ceramic,  the  value  of  the  HEL,  the 
longitudinal  modulus  and  the  underlying  hydrostat  can  be  different.  Loss  of  strength  due  to 
damage  in  the  ceramic  will  collapse  the  Hugoniot  toward  the  hydrostat  of  the  fully  dense 
ceramic.  With  incomplete  crushing  out  of  porosity,  the  ceramic  will  have  an  unload  cold 
compressive  hydrostat  determined  by  the  specific  volume  on  unload  that  is  intermediate 
between  theoretical  density  and  the  initial  density.  Energy  expended  to  crush  out  the  porosity 
will  become  internal  energy  and  also  affect  the  shape  of  the  hydrostat.  Thus,  if  the  porosity 
is  continually  changing  due  to  the  stress  level  then  the  hydrostat  would  also  be  continually 
changing.  This  effect  in  ceramics  can  totally  confuse  the  interpretation  of  experimental 
results  and  is  difficult  to  model  accurately. 

The  effects  that  the  RDA  ceramic  model  has  on  the  particle  velocity  history  are  best 
understood  in  the  context  of  the  idealized  load/unload  wave  profile  shown  in  Figure  3.  This 
figure  is  a  schematic  of  an  experimental  planar  impact  geometry.  Position  0.0  is  the  flyer 
plate/target  plate  interface,  the  left  boundary  is  the  flyer  plate  and  the  right  boundary  is  the 
target  ceramic/backplate  interface.  The  planar  impact  produces  a  compressive  wave  of 
uniaxial  strain  which  propagates  through  both  plates.  If  the  flyer  and  target  plate  are  the 
same  material  (as  shown  in  the  figure)  the  wave  is  symmetric  about  the  impact  plane.  At 
high  impact  velocities,  the  shock  wave  will  be  characterized  by  a  two  wave  structure.  Until 
the  stress  in  the  material  reaches  its  elastic  limit  (in  uniaxial  geometry,  the  HEL)  the  wave 
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leading  the  shock  is  elastic  and  travels  at  the  elastic  wave  speed.  This  speed  is  strictly  a 
function  of  the  elastic  constants,  p0,  C0,  and  v  and  can  be  determined  pnor  to  the  experiment. 
The  time  at  which  this  wave  arrives  at  the  backplate  is  used  for  calibrating  of  the  elastic 
speed.  Above  the  elastic  limit  the  shock  wave  is  characterized  by  the  second  or  ’plastic’ 
wave.  The  arrival  time  of  this  wave  is  dependent  upon  the  strength  of  the  material  and  the 
hydrostatic  pressure.  The  shape  or  profile  of  the  plastic  wave  front  is  also  determined  by  the 
strength  and  pressure. 

The  two  symmetric  elastic  and  plastic  shock  waves  travel  in  opposite  directions 
eventually  reaching  an  interface  and  unloading.  The  unload  is  complete  (zero  stress)  if  the 
interface  is  a  free  surface.  Conversely,  if  the  interface  is  a  perfect  impedance  r  atch  there  is 
total  transmission  of  the  wave  and  no  unload  wave.  Therefore  any  impedance  at  the 
boundary  that  is  comparable  to  the  ceramic  impedance  will  affect  the  magnitude  of  the 
unloading  wave.  The  shape  and  speed  of  the  unloading  wave  are  equally  important  for 
modeling  purposes.  The  wave  first  unloads  elastically.  The  magnitude  of  this  unload  is 
directly  related  to  the  strength  of  the  material  and  is  equal  to  4/3  Y,  where  Y  is  the  strength. 
Therefore  any  difference  between  the  strength  used  to  predict  the  loading  wave  and  the 
unloading  wave  can  be  ascribed  to  further  ’softening’  such  as  damage  that  results  after  load. 
The  second  wave  is  the  plastic  unloading  wave  and  like  the  plastic  loading  wave  is  dependent 
on  the  strength  and  pressure. 

The  unloading  wave  travels  in  the  opposite  direction  as  the  loading  wave.  Because  the 
flyer  plate  is  thinner  than  the  target  plate  in  the  experiments  modeled,  the  two  unloading 
waves  intersect  in  the  target  plate.  With  totally  free  surfaces,  unlimited  strength,  and  in  a 
lossless  medium  the  magnitude  of  the  stress  of  these  intersecting  waves  would  be  negative 
and  twice  the  shock  stress.  Normally  the  material  fails  (spalls)  long  before  this  stress  is 
reached. 

For  the  initial  calculations  for  predicting  the  measured  particle  velocity  history,  a 
functional  form  for  the  hydrostat  was  assumed,  and  the  ceramic  was  treated  as  behaving 
perfectly  plastic.  Then  the  only  active  parameter  in  the  Johnson-Cook  flow  stress  model  is 
Cj,  the  initial  yield  stress.  To  chose  an  initial  yield,  the  well  known  relationship 
between  dynamic  yield  and  HEL  was  used: 
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where 


~  °HEL  d  ■  2v)/(l  -  v) 


(6) 


Orel  =  Hugoniot  elastic  limit 
v  =  Poisson’s  ratio 

The  HEL  value  was  calculated  from  the  particle  velocity  history  at  the  ceramic/backplate 
interface  by  the  expression 

Ohel  -  ^  Um  P) 

2 

where  Zc  and  Zj  are  the  appropriate  shock  impedance  for  the  ceramic  and  backplate, 
respectively,  and  u,,,  is  the  observed  particle  velocity  amplitude  selected  from  each  profile 
which  represents  the  transition  from  elastic  to  nonelastic  behavior. 

As  will  be  shown,  for  the  materials  modeled,  this  approach  accurately  matched  the 
break  point  in  the  velocity  history  at  the  HEL  and  the  peak  stress  level.  This  approach  is 
normally  deficient  at  predicting  1)  the  profile  of  the  plastic  shock  front  on  load,  2)  the 
magnitude  of  the  elastic  unload  if  the  material  has  softened  as  a  result  of  compression,  and  3) 
the  profile  of  the  plastic  unload.  It  then  becomes  necessary  to  make  adjustments  within  the 
parameter  space  of  the  yield  function.  In  this  analysis,  it  was  assumed  that  work  hardening, 
strain  rate  effects,  and  damage  softening  were  the  relevant  parameters  to  adjust.  The  effects 
of  thermal  softening,  and  confining  pressure,  which  are  also  in  the  model,  have  not  been 
included  in  the  present  analysis.  Work  hardening  is  described  with  a  power  law  function  (see 
equation  3),  noting  that  the  exponent  characterizes  the  sensitivity  of  strength  to  plastic  strain. 
Variation  in  the  work  hardening  coefficients  has  the  effect  of  varying  the  apparent  HEL  and 
increases  in  work  hardening  stress  will  increase  the  HEL.  Consequently,  to  preserve  the  HEL 
indicated  by  the  experimental  results,  the  static  yield  stress  must  be  reduced.  The  strain  rate 
term,  C3,  modifies  the  net  strength  determined  by  the  initial  yield  and  work  hardening  terms. 
Strain  rate  effects  can  also  increase  the  HEL  and  magnitude  of  the  elastic  unload. 

Application  of  the  RDA  compressive  damage  model  acts  to  reduce  the  strength  of  the 
material.  Damage  softening  has  been  utilized  to  a)  reduce  the  arrival  time  of  the  plastic 
wave,  b)  modify  the  profile  of  the  plastic  wave,  and  c)  reduce  the  magnitude  of  the  elastic 
unload. 
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As  both  the  Johnson-Cook  model  and  the  RDA  compressive  damage  model  are 
phenomenological,  confidence  of  their  accuracy  is  rooted  in  their  ability  to  predict  results 
over  impact  velocity  and  spatial  geometry  regimes  of  interest. 
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4.  THE  PC-WONDY  COMPUTER  CODE: 


Since  the  requirements  for  computational  speed  and  memory  are  gready  reduced  for 
problems  involving  only  one  dimension  it  was  decided  early  on  in  the  program  to  identify  a 
hydro  code  that  would  lend  itself  to  easy  porting  to  a  PC,  allow  the  implementation  of  the 
RDA  ceramic  model,  and  be  suitable  for  the  types  of  problems  of  interest,  such  as  plate 
impact  experiments.  The  WONDY-V  hydrocode  is  a  one-dimensional  finite-difference  wave 
propagation  code  (Kipp  and  Lawrence,  1981)  and  was  chosen  as  the  vehicle  for  the  modeling 
effort. 

WONDY  is  comprised  of  about  four  thousand  lines  of  FORTRAN  making  its  code 
core  size  small,  even  for  a  personal  computer.  WONDY  maintains  all  zonal  quantities  for  all 
zones  in  one,  singly  subscripted  array  and  uses  indirect  memory  access  to  store  and  retrieve 
zonal  information.  This  approach  has  a  number  of  advantages  over  codes  that  maintain  zonal 
quantities  in  separate  arrays.  First,  it  allows  for  ease  of  splitting  and  recombining  zones 
should  rezoning  be  a  requirement  of  a  particular  problem  (e.g.  tracking  a  very  steep  shock). 
Second,  the  size  of  output  "tapes"  can  be  kept  small  since  the  data  requested  for  output  can 
be  identified  by  an  index  (as  opposed  to  placing  all  variables  out  on  a  tape).  More 
importantly,  it  provides  for  easy  implementation  of  material  models  that  have  an  arbitrary 
number  of  state  variables.  WONDY  also  has  available  an  extensive  library  of  equation-of- 
state  (EOS)  routines  useful  for  characterizing  a  number  of  solids,  gases,  and  explosives.  It  is 
felt  that  WONDY’s  treatment  of  the  hydrodynamics  and  thermodynamics  of  impact  is 
superior  to  other  available  codes.  In  particular  the  advantage  of  using  local  densities  and 
sound  speeds  can  not  be  over  emphasized  when  trying  to  predict  accurate  deviatoric  stress 
levels  and  shock  propagation  velocities. 

The  environment  of  the  PC  provides  the  user  with  an  extensive  array  of  software 
utilities,  graphics  libraries,  etc.  that  can  boost  both  the  productivity  of  code  development  and 
quality  of  graphic  output.  Making  WONDY  part  of  this  environment  has  significantly  aided 
the  ability  to  understand  the  response  of  ceramics  to  dynamic  loading. 
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5.  SANDIA  PLATE  IMPACT  EXPERIMENTS: 


The  experiments  selected  for  simulation  were  the  silicon  carbide  and  boron  carbide 
plate  impact  experiments  recently  reported  by  Kipp  and  Grady  (1989).  Note  that  this  report 
also  contains  results  for  titanium  diboride  and  zirconium  dioxide.  A  typical  test  setup  is 
shown  in  Figure  4.  These  impact  experiments  were  symmetric,  meaning  that  the  flyer  and 
target  plate  are  the  same  material.  Symmetric  impact  experiments  have  the  implicit  benefit  of 
eliminating  any  uncertainties  associated  with  the  modeling  of  a  second  material.  However, 
two  other  materials  necessary  for  the  experiment  appear  in  the  figure  and  must  be  accounted 
for  when  modeling  the  entire  experiment.  The  first  was  a  foam  backing  material  used  to  push 
the  flyer  plate.  Depending  on  the  anticipated  flyer  velocity,  foam  with  two  different  densities 
was  used.  The  lower  density  foam  was  0.32  g/cc  and  the  higher  density  foam  was  0.64  g/cc 
and  was  therefore  on  the  order  of  one-quarter  the  density  of  the  boron  carbide.  On  the  back 
side  of  the  target  plate  was  mounted  a  lithium  fluoride  crystal  window.  Laser  velocity 
interferometry  techniques  were  used  to  resolve  the  longitudinal  motion  at  the  ceramic-lithium 
fluoride  interface.  Both  the  density  and  sound  speed  of  the  lithium  fluoride  are  comparable  to 
those  of  the  ceramics  modeled  and  therefore  had  a  significant  effect  on  the  magnitude  of  the 
transmitted  wave  through  the  ceramic-lithium  fluoride  interface.  The  interface  particle  motion 
versus  time  was  the  principal  diagnostic  of  the  experiment  and  the  measure  by  which  to  judge 
the  accuracy  of  the  model. 

Plate  dimensions  were  chosen  with  two  criteria:  1)  the  release  wave  coming  from  the 
foam-ceramic  interface  would  intersect  the  ceramic-lithium  fluoride  release  wave  in  the  target 
plate  and  2)  release  waves  from  the  lateral  boundaries  of  the  plates  would  not  reach  the  axis 
of  measurement  until  after  times  of  interest.  In  addition,  measurements  were  made  on  each 
ceramic  with  two  flyer  plate  velocities,  designated  a  low  and  high  velocity.  In  general  the 
low  velocity  case  was  about  1500  m/s  and  the  high  velocity  case  was  about  2100  m/s.  The 
particle  velocity  versus  time  results  for  the  SiC  and  the  B4C  are  shown  in  Figures  5  and  6, 
respectively. 
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6.  MODEL  RESULTS: 


6.1  Ceramic  Properties: 

The  two  ceramic  materials  that  were  modeled  are  SiC  and  B4C.  Table  1  shows  the 
elastic  properties  for  SiC  and  B4C  (Kipp  and  Grady,  1989)  that  were  used  in  the  RDA 
ceramic  model.  Figure  7  shows  the  Hugoniot  for  SiC  using  LANL  data  (Marsh,  ed.,  1980). 
The  solid  line  is  a  fit  to  the  data  for  t]  greater  that  0.15  assuming  a  linear  shock  versus 
particle  velocity  relationship  and  represents  an  assumed  shape  of  the  underlying  hydrostat 
from  equation  5.  The  dashed  lines  represent  a  classical  elastic  perfectly-plastic  behavior.  SiC 
appears  to  closely  resemble  this  type  behavior  from  the  apparent  HEL  at  120  kbars  to 
approximately  420  kbars.  Note  the  apparent  loss  of  strength  for  1]  values  of  0.14  to  0.16  and 
the  collapse  to  the  hydrostat. 

The  Hugoniot  data  for  B4C  is  much  more  complex.  Shown  in  Figure  8  is  Hugoniot 
data  (Gust  et  al.,  1971).  Again  the  solid  line  is  a  fit  to  the  data  for  q  greater  than  0.15  as 
was  done  for  the  SiC.  The  dashed  lines  represent  the  classical  elastic  perfectly-plastic 
behavior.  Substantial  scatter  in  the  data  around  the  apparent  HEL  at  150  kbars  is  observed. 
The  loading  curve  to  the  two  impact  velocities  modeled  (shown  with  the  vertical  dashed  lines) 
is  far  from  plastic  behavior  and  indicates  loss  of  strength  and  possible  changes  in  the 
hydrostat. 


TABLE  1.  CERAMIC  PROPERTIES 


Material 

po  , 

CL 

cs 

C0 

V 

°HEL 

kg/m 

km/s 

km/s 

km/s 

GPa 

SiC 

3177 

12.06 

7.67 

8.19 

0.16 

15.0 

b4c 

2516 

14.04 

8.9 

9.57 

0.164 

14.4 
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6.2  Model  Results  for  SiC: 


The  low  and  high  velocity  impact  particle  velocity  experimental  results  for  SiC  are 
shown  in  Figure  5.  Of  the  two  materials  modeled,  the  silicon  carbide  more  closely  resembles 
a  classical  elastic-plastic  behavior  as  illustrated  in  Figure  1.  The  elastic  and  plastic  shock 
wave  structure  and  the  elastic  and  plastic  unload  structure  are  both  readily  apparent  An 
additional  wave,  most  apparent  in  the  high  velocity  case  on  load  and  unload  was  caused  by 
the  wave  reflected  at  the  lithium  fluoride  window.  A  comparison  of  the  measured  and 
predicted  results  using  a  perfectly  plastic  model  for  the  low  and  high  velocity  measurements, 
respectively,  is  shown  in  Figure  9(a  &  b).  This  simple  model  does  a  remarkably  good  job  of 
predicting  the  measured  waveforms,  particularly  in  the  high  velocity  impact  case. 

Both  simulations  lack  the  necessary  yield  function  to  fill  the  ’notched’  region  on  load. 
The  perfectly  plastic  model  for  the  low  velocity  case  over-predicts  the  level  of  the  plastic 
unload.  The  initial  hypothesis  was  that  the  material  spalled  and  with  appropriate  adjustment 
of  the  tensile  damage  parameters  could  be  adequately  modeled. 

Attention  was  first  focused  on  modeling  the  obvious  features  of  load  and  elastic 
unload  by  variation  of  the  strain  rate,  strain  hardening  and  damage  coefficients.  Figure 
10(a-d)  shows  a  progression  of  computations  varying  only  the  strain  rate  coefficient,  C3  for 
the  high  velocity  case.  Choosing  a  value  of  .02  for  C3  appears  to  yield  an  excellent  match  to 
the  observed  velocity  versus  time  profile  for  features  related  to  the  load  and  elastic  unload. 
But  there  are  subtle  features  that  cannot  be  matched  with  strain  rate  alone.  This  choice  for 
C3  removed  too  much  of  the  comer  at  the  loading  peak.  Also,  the  magnitude  of  the  elastic 
unload  is  too  large.  The  unload  wave  arrives  too  early,  suggesting  that  a  dispersive 
mechanism  has  occurred  in  the  plate.  Finally,  even  though  the  plastic  unload  level  appears 
correct  there  was  an  anomalous  wave  that  arrives  at  1.7  p.s. 

Figure  ll(a-d)  shows  the  collective  effort  of  many  calculations  resulting  in  choices  for 
all  parameters  that  in  an  overall  sense  yield  the  best  fit  to  the  load  and  the  elastic  unload 
features  of  the  high  velocity  data.  In  order,  the  progression  is  from  perfectly  plastic  with  an 
initial  yield  value  of  Cj  =  80  kbar,  to  the  addition  of  75  -  80  kbar  of  work  hardening  using 
values  C2  =  100,  n  =  0.1,  to  the  addition  of  6  -  8  percent  compressive  damage  using  a  value 
a  =  0.7  x  (lO)"6  (kg/J),  to  the  addition  of  a  strain  rate  coefficient  of  0.008.  Both  work 
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hardening  and  strain  rate  have  the  effect  of  filling  in  the  ’notched*  region  seen  in  the  perfectly 
plastic  case.  However,  strain  rate  tends  to  overfill  this  region  near  the  HEL  and  under  fill  the 
region  at  the  peak.  The  work  hardening  function  fills  this  region  more  uniformly  and  overall 
allows  for  more  accurate  matching  to  the  profile  of  the  plastic  wave. 

In  all  cases  the  unloading  wave  arrives  too  soon.  In  order  to  match  this  arrival  it  was 
necessary  to  modify  the  isentrope  during  unload.  The  conventional  approach  to  using  a  Mie- 
Gruniesen  EOS  in  a  hydro  code  is  to  unload  along  the  load  isentrope.  It  should  be  noted  that 
the  velocity  histories  reported  by  Kipp  and  Grady  indicated  significant  dispersion  observed  in 
the  unload  for  boron  carbide,  titanium  diboride,  and  zirconium  dioxide  materials  with  SiC 
indicating  somewhat  less  dispersion.  If  the  unload  isentrope  was  taken  along  a  path  to  a 
lower  density  than  the  original  density  at  zero  pressure  then  a  suitable  match  with  the  SiC 
could  be  made.  PC-WONDY  was  modified  to  allow  specification  of  an  arbitrary  pressure 
offset  at  zero  strain  during  unload.  Without  suitable  theory  to  support  such  an  effect  caution 
has  been  exercised  before  using  this  option.  However,  the  model  results  are  improved  when  a 
30  kbar  offset  at  zero  strain  was  used  as  shown  in  Figure  12  which  compares  the  model  with 
the  measurement. 

Shown  in  Figure  13  is  the  comparison  between  the  model  prediction  and  the  low 
velocity  measurement  using  the  same  parameters  as  in  the  high  velocity  case.  These  results 
compare  favorably  with  the  experimental  data  for  the  first  and  second  shock  wave  on  load, 
and  the  arrival  time  of  the  elastic  unload  wave.  However,  the  calculation  predicts  only  elastic 
unload  whereas  the  measurement  shows  a  plastic,  highly  dispersive  wave  at  1.35  ps. 

The  arrival  of  the  additional  wave  at  1.7  ps  and  the  modest  dispersive  unload  wave  in 
the  high  velocity  case  and  the  highly  dispersive  unload  wave  in  the  low  velocity  case  were  all 
originally  thought  to  be  a  spall  phenomenon  that  could  be  modeled  by  suitable  adjustment  of 
the  tensile  damage  model  parameters.  Figure  14  is  used  to  show  the  characteristic  particle 
velocity  history  that  occurs  when  a  material  spalls.  The  abrupt  change  in  slope  seen  in  the 
figure  is  not  observed  in  either  the  low  or  high  speed  impact  particle  velocity  measurements 
suggesting  the  material  either  has  no  spall  strength  or  does  not  spall  and  some  other  event 
was  responsible  for  the  wave  structure.  Further,  the  model  predicted  no  spall,  regardless  of 
the  magnitude  of  the  tensile  model  damage  sensitivity  parameters.  Closer  examination  of  the 
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calculation  revealed  that  the  mean  stress  (pressure)  remained  positive  (compressive) 
throughout  the  entire  calculation.  Figure  15  shows  the  pressure  distribution  when  the  pressure 
reaches  its  absolute  minimum  in  the  target  plate.  Since  the  tensile  model  is  activated  by  a 
negative  pressure,  no  tensile  damage  was  calculated.  The  ceramic  maintained  this 
compressive  state  for  three  reasons:  1)  the  impedance  of  lithium  fluoride  was  sufficiently 
close  to  the  ceramic  target  to  allow  transmission  of  a  major  fraction  of  the  pressure 
component  of  the  shock  wave  (using  the  computed  peak  particle  velocities  and  equations  (2- 
7)  to  (2-9),  the  impedances  of  the  silicon  carbide  and  lithium  fluoride  are  4.3  x  107  and  1.62 
x  107  rayls,  respectively),  2)  30  percent  of  the  pressure  was  transmitted  into  the  foam  backing 
material,  and  3)  energy  was  deposited  in  the  plates  that  acts  to  reduce  the  strength  of  the 
ceramic  either  through  damage  (micro-cracking)  or  thermal  melting.  For  the  modeled  impact 
velocities,  removal  of  the  lithium  fluoride  window  is  sufficient  to  produce  negative  pressure 
in  the  target  and  cause  spall.  For  comparison,  Figure  16  shows  the  pressure  distribution  for 
the  high  velocity  case  when  the  lithium  fluoride  window  has  been  removed  from  the 
calculation.  Note  the  signature  of  a  spall  zone  at  a  position  6  to  7  mm  from  the  flyer/target 
plate  interface. 

If  the  material  indeed  does  not  spall  then  another  explanation  is  required  for  the  wave 
structure  near  1.7  |is  after  plate  impact.  Further  investigation  with  variation  of  the  boundary 
conditions  at  interfaces  including  removing  the  foam  material  from  the  calculation  revealed 
that  the  arrival  time  of  the  anomalous  wave  was  coincident  with  the  arrival  of  a  wave  caused 
by  an  event  occurring  at  the  boundary  of  the  foam-ceramic  interface  immediately  after 
unload.  This  result  lead  to  the  conclusion  that  the  observed  wave  form  is  not  a  spall  signal 
and  supports  the  model  result  that  the  target  ceramic  did  not  spall.  Therefore,  some 
modification  to  the  damage  model  is  necessary  to  produce  the  observed  events  associated  with 
unload  in  the  flyer. 

The  comparison  between  the  model  prediction  and  measured  particle  velocity  for  the 
low  impact  SiC  shows  that  the  model  over-predicts  plastic  unload  whereas  the  measurement 
shows  a  highly  dispersive  plastic  unload  wave.  Like  the  high  speed  case  the  low  speed  case 
does  not  show  the  characteristic  spall  signal  as  illustrated  in  Figure  14.  The  unload  level  is 
directly  related  to  the  amount  energy  released  during  unload  i.e.  energy  transmitted  into  the 
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foam  or  dissipated  in  the  flyer  or  target  will  tend  to  reduce  the  amount  of  unload  increasing 
the  particle  velocity  level.  Since  energy  lost  to  micro-cracking  is  known  to  reduce  the 
magnitude  of  the  elastic  unload,  the  compressive  damage  model  can  be  used  as  a  device  to 
reduce  the  strength  of  the  material,  beyond  what  was  necessary  to  match  load,  to  see  what 
effect  this  would  have  on  the  unload  particle  velocity  history.  Figure  17(a  -  d)  shows  the 
progression  of  results  when  a  was  systematically  increased.  Figure  18  compares  the 
experimental  measurements  with  the  model  results  using  an  a  =  0.32  x  (10)'5  (kg/J)  which 
best  matches  the  elastic  unload.  This  value  of  damage  sensitivity  results  in  damage  in  the 
target  plate  of  40  percent.  Obviously  increasing  a  ruins  the  predicted  performance  for  load. 
Since  both  the  load  and  unload  waves  cannot  be  matched  with  the  same  value  of  a  and  since 
the  load  can  be  correctly  predicted,  then  damage  or  some  other  strength  reducing  mechanism 
was  taking  place  after  load  while  the  material  was  at  peak  stress  or  even  as  it  unloads. 
Candidate  theories  are  fatigue,  melting  at  grain  boundaries  of  the  ceramic  and  conversion  of 
internal  energy  to  micro-crack  growth  and  coalesce.  Since  the  model  predicts  the  high  speed 
impact  adequately  without  introducing  additional  damage,  the  additional  damage  may  be 
inversely  related  to  impact  speed.  This  suggests  a  confuting  pressure  effect,  or  other  crack 
growth  mitigating  effect  such  as  ratio  between  pressure  and  deviatoric  stress. 

Silicon  carbide  has  been  successfully  modeled  with  the  RDA  ceramic  failure  model 
with  the  Johnson-Cook  flow  stress  model  modified  by  damage  softening.  Within  the  context 
of  this  model  SiC  exhibits  about  75-80  kbar  of  work  hardening.  A  compressive  damage 
sensitivity  of  0.7  x  (10)-6  (kg/J)  produced  6-8  percent  and  12  -  15  percent  strength  reduction 
on  load  in  the  low  and  high  velocity  impact  cases,  respectively.  Some  strain  rate  dependence 
was  required  to  smooth  certain  portions  of  the  velocity  history.  It  was  also  necessary  to 
prescribe  a  separate  unloading  isentrope  with  a  pressure  offset  of  30  kbar  to  match  the  arrival 
time  of  the  unload  wave  in  the  high  speed  case. 

Model  predictions  made  a  satisfactory  match  with  experiment  for  load,  but  over¬ 
predicted  the  magnitude  of  unload  for  the  low  speed  impact  case.  After  loading,  increasing 
the  calculated  damage  from  the  predicted  level  of  6  -  8  percent  to  40  percent  while  the 
material  was  at  peak  stress  was  necessary  to  match  the  unload  particle  velocity  history. 


16 


6.3  Model  Results  for  I^C: 

The  approach  to  modeling  B4C  begins  with  an  examination  of  the  available  Hugoniot 
data.  Figure  19  is  stress  versus  T|  for  B4C  from  planar  impact  experiments  (Gust,  et  al, 

1971).  The  solid  line  is  a  fit  to  the  data  for  r\  greater  than  0.15  and  is  the  shape  of  the 
hydrostat  assuming  a  linear  relationship  between  the  shock  velocity  and  the  particle  velocity 
(equations  6  and  7).  The  observed  Atj  at  constant  stress  for  q  =  0.055  to  .065  is  indicative 
of  possible  loss  of  porosity  or  phase  change  in  the  ceramic.  The  shape  of  the  underlying 
hydrostat  was  therefore  unknown  and  the  load  and  unload  paths  would  be  different. 

The  low  and  high  velocity  impact  particle  velocity  experimental  results  are  shown  in 
Figure  6.  Unlike  the  SiC,  B4C  exhibits  a  number  of  differences  distinguishing  its  behavior 
from  the  ideal  elastic-plastic  behavior  shown  in  Figure  3.  First,  the  HEL  is  marked  by 
uncharacteristic  oscillations.  Kipp  and  Grady  account  for  this  behavior  as  a  "heterogeneous 
failure  or  faulting  mechanism  during  compressive  loading  coupled  with  the  small  laser  spot 
size  associated  with  VISAR  measurements."  No  better  explanation  is  apparent,  nor  does  the 
model  attempt  to  predict  it.  Instead,  values  for  the  yield  parameters  values  have  been  chosen 
that  best  fit  the  mean  of  the  oscillatory  behavior.  Second,  completely  absent  in  the  particle 
velocity  measurement  of  the  low  speed  impact  case  was  any  evidence  of  elastic  unload.  In 
the  classical  case  there  is  always  a  linear  non-dispersive  unloading  path.  This  absence  of 
elastic  unload  suggests  a  total  loss  of  strength.  Third,  the  entire  unloading  path  was 
dispersive  and  cannot  be  accounted  for  by  the  hydrostatic  model  used  during  load.  This 
suggests  that  the  material  has  substantially  changed  in  character,  perhaps  undergoing  phase  or 
phases  of  transformation  and/or  as  with  SiC,  the  material  has  additional  damage  or  other 
strength  reducing  mechanism  that  was  taking  place  after  load  while  the  material  was  at  peak 
stress  or  even  as  it  unloads. 

In  the  low  impact  velocity  measurement,  a  dip  was  observed  immediately  following 
the  arrival  of  the  plastic  wave  when  one  normally  expects  to  see  a  plateau.  Although  it  is  not 
unusual  behavior,  the  markedly  higher  sound  speed  of  the  B4C  created  some  modeling 
diagnostic  difficulty  in  the  context  of  the  plate  thicknesses  and  impact  velocities  used  in  the 
experiments.  The  combination  of  a  slow  plastic  loading  wave  (caused  by  considerable 
damage)  and  fast  elastic  unloading  wave  allowed  the  unloading  wave  to  overtake  the  loading 
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wave  at  approximately  the  VISAR  location.  This  event  was  further  complicated  by  the  wave 
reflecting  at  the  ceramic -lithium  fluoride  interface. 

Shown  in  Figure  20  is  a  comparison  of  the  measured  and  predicted  results  using  an 
elastic  perfectly-plastic  model  for  the  low  velocity  impact  case.  The  model  can  adequately 
predict  the  HEL,  but  the  similarity  ends  there.  In  the  perfectly  plastic  model  the  plastic  shock 
arrival  time  was  too  early.  The  measurements  clearly  indicate  a  transition  from  an  elastic  to  a 
perfectly  plastic  wave.  The  addition  of  damage  with  the  compressive  damage  model  was  then 
used  to  delay  the  arrival  of  the  plastic  shock.  Figure  21(a-d)  shows  the  predicted  waveforms 
after  making  adjustments  of  the  compressive  damage  sensitivity  parameter,  a  (see  equation 
1  -  20)  for  the  low  velocity  impact  case.  Note  in  the  figures  that  as  a  increases,  the  arrival 
time  of  the  plastic  shock  was  delayed.  The  elastic  unloading  wave,  whose  arrival  time  will 
not  vary  (so  long  as  the  material  retains  some  elasticity),  begins  to  catch  up  to  the  plastic 
load,  eventually  overtaking  it  in  Figure  21(d). 

Figures  22  and  23  show  the  comparison  between  the  model  predictions  and  the 
measured  particle  velocity  histories  for  both  the  low  and  high  velocity  cases.  These  results 
were  obtained  using  a  damage  parameter  of  a  =  3.4  x  (10)"5  kg/J  and  produces  a  good  fit  to 
the  observed  loading  waveforms.  The  level  of  damage  computed  as  a  result  of  the 
compressive  shock  was  between  45  -  50  percent,  and  65  -  70  percent  for  the  low  and  high 
velocity  cases,  respectively.  As  mentioned  earlier,  the  measurements  for  the  low  and  high 
velocity  impact  experiments  indicate  highly  dispersive  unloading  waves.  The  low  speed 
impact  case  shows  no  evidence  of  any  elastic  strength  after  being  compressed.  The  small 
linear  segment  at  the  beginning  of  unload  in  the  high  speed  impact  case  was  believed  to  be 
evidence  of  the  residual  strength  of  the  ceramic.  Using  the  length  of  this  segment,  the 
measured  particle  velocity  and  equations  (2  -  7)  to  (2  -  9),  implies  a  residual  strength  of  15.4 
kbar,  implying  damage  of  87  percent  in  the  ceramic.  The  linear  segments  at  the  beginning  of 
unload  in  the  low  and  high  speed  plate  impact  model  predictions  indicates  that  the  material  in 
the  two  experiments  should  have  more  strength  than  was  observed. 

Boron  carbide  has  been  successfully  modeled  by  the  RDA  ceramic  failure  model  with 
the  Johnson-Cook  flow  stress  model  modified  by  damage  softening.  B4C  exhibits  little  or  no 
work  hardening,  and  no  significant  strain-rate  dependence.  Using  a  value  of  a  =  3.4  x  (10)'5 
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kg/J  provides  a  good  fit  to  the  measured  loading  waves  for  both  low  and  high  velocity  cases. 
The  level  of  damage  computed  as  a  result  of  the  compressive  shock  was  between 
45  -  50  percent,  and  65  -  70  percent  for  the  low  and  high  velocity  cases,  respectively. 
Dispersive  unloading  with  little  evidence  of  an  elastic  unload  are  required  to  predict  the 
experimental  data  suggesting  either  damage  reduced  modulus  or  non-zero  pressure  on  return 
to  original  density.  Removing  the  strength  from  the  ceramic  in  the  model  predictions  to 
match  the  unload  profile  affects  the  good  match  with  the  load  profile  and  leads  to  the  same 
conclusion  that  was  reached  with  the  SiC.  Damage  or  some  other  strength  reducing 
mechanism  was  at  work  after  load,  decreasing  the  strength  of  the  ceramic  more  than  was 
necessary  to  accurately  predict  load.  A  further  similarity  between  the  two  materials  was  that 
the  damage  appears  to  be  more  severe  at  the  low  impact  velocity  than  at  the  high  impact 
velocity. 
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7.  SUMMARY 


The  RDA  ceramic  failure  model  was  used  to  model  the  behavior  of  SiC  and  B4C 
during  plate  impact  loading.  Experiments  modeled  were  those  of  Sandia  National  Laboratory 
(Kipp  and  Grady,  1989).  The  simulations  were  performed  with  a  PC  based  version  of  the 
Sandia  WONDY-V  hydrocode  (Kipp  and  Lawrence,  1981).  WONDY-V  was  converted  to  a 
386  PC,  and  then  the  RDA  ceramic  failure  model  was  incorporated  into  it. 

The  model  characterized  the  two  ceramic  materials  with  distinctly  different  work 
hardening,  strain  rate,  and  damage  sensitivity.  Within  the  context  of  the  model  SiC  exhibited 
75  -  80  kbar  work  hardening  and  some  strain  rate  sensitivity.  A  compressive  damage 
sensitivity  of  0.7  x  (10)"6  (kg/J)  produced  6-8  percent  and  12  -  15  percent  strength  reduction 
on  load  in  the  low  and  high  velocity  impact  cases,  respectively.  Some  strain  rate  dependence 
was  required  to  smooth  certain  portions  of  the  velocity  history. 

Model  predictions  for  SiC  made  a  satisfactory  match  with  experiment  for  load,  but 
over-predicted  the  magnitude  of  unload  for  the  low  speed  impact  case.  After  loading, 
increasing  the  calculated  damage  from  the  predicted  level  of  6  -  8  percent  to  40  percent  while 
the  material  was  at  peak  stress  was  necessary  to  match  the  unload  particle  velocity  history. 

In  contrast  to  SiC,  B4C  exhibits  little  or  no  work  hardening,  and  no  significant  strain- 
rate  dependence.  Using  a  value  of  a  =  3.4  x  (10)'5  kg/J  provides  a  good  fit  to  the  measured 
loading  waves  for  both  low  and  high  velocity  cases.  The  level  of  damage  computed  as  a 
result  of  the  compressive  shock  was  between  45  -  50  percent,  and  65  -  70  percent  for  the  low 
and  high  velocity  cases,  respectively.  Dispersive  unloading  with  little  evidence  of  an  elastic 
unload  is  required  to  predict  the  experimental  data  suggesting  either  damage  reduced  modulus 
or  non-zero  pressure  on  return  to  original  density. 

For  both  ceramics,  particle  velocity  histories  indicate  that,  after  loading,  further 
damage  occurs  and  is  inversely  related  to  impact  velocity.  When  damage  is  allowed  to 
increase  beyond  the  level  required  to  adequately  match  the  loading  profile,  a  suitable  match 
with  the  unloading  profile  can  be  made.  This  additional  damage  which  cannot  be  predicted 
with  the  current  model  suggests  other  strength  reducing  mechanisms  are  at  work  while  the 
ceramic  was  at  peak  stress  or  at  unload.  Candidate  theories  for  this  behavior  are  fatigue 
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and/or  thermal  melting  at  the  grain  boundaries  of  the  ceramic.  Lack  of  additional  damage 
evolution  in  the  high  speed  impact  case  suggests  its  cause  may  be  mitigated  by  confining 
pressure  or  ratio  of  pressure  to  deviatoric  stress  effects. 

Positive  (compressive)  mean  stress  levels  were  maintained  throughout  the  test 
specimen  during  the  experimental  measurement  as  a  result  of  the  boundary  conditions 
provided  by  the  lithium  fluoride  window  of  the  VISAR  instrumentation  and  the  flyer  plate 
foam  backing  material.  Energy  absorption  as  a  result  of  compressive  damage  and/or  melting 
at  the  grain  boundaries  also  contributes  to  the  maintaining  of  positive  mean  stress  levels.  The 
predicted  axial  stress  level  in  the  high  impact  velocity  case  for  SiC  was  approximately  25 
kbar  which  is  well  above  the  3-5  kbar  range  normally  associated  with  the  spall  strength  of 
low  porosity,  high  strength  ceramics,  however,  no  spall  was  evident  in  the  experimentally 
measure  signals  and  the  model  predicts  no  spall  to  occur.  Therefore,  these  results  support  the 
hypothesis  that  for  tensile  failure  in  these  ceramics,  the  mean  stress  must  be  negative  (tensile). 
This  result,  if  true,  would  have  a  significant  impact  on  the  design  of  armor.  A  major  design 
issue  would  be  whether  or  not  the  ceramic  is  confined,  not  the  level  of  confining  pressure. 


THE  VIEWS,  OPINIONS,  AND/OR  FINDINGS  CONTAINED  IN  THIS  REPORT  ARE 
THOSE  OF  THE  AUTHORS  AND  SHOULD  NOT  BE  CONSTRUED  AS  AN  OFFICIAL 
DEPARTMENT  OF  THE  ARMY  POSITION,  POLICY,  OR  DECISION,  UNLESS  SO 
DESIGNATED  BY  OTHER  DOCUMENTATION. 
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STRESS 


FIGURE  1: 


STRAIN  (%) 


Elastic  perfectly-plastic  stress-strain  path  showing  relationship  to  the  hydrostat 
(heavy  line).  Dashed  line  illustrates  response  with  softening  such  as  damage. 
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FIGURE  2:  Schematic  of  Hugoniot  for  a  ceramic  material  with  porosity  indicating  that  at  each 
initial  specific  volume,  the  material  will  have  a  different  HEL  and  potentially 
different  longitudinal  modulus.  Loss  of  strength  due  to  damage  in  the  material 
collapses  the  Hugoniot  toward  the  hydrostat  of  the  fully  dense  material. 
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AFTER  LOAD 


AFTER  UNLOAD 


FIGURE  3:  Idealized  load/unload  wave  profiles  for  planar  impact  experiments. 
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SUPPORT  RING 


FIGURE  4:  Plate  impact  experimental  configuration  of  Sandia  National  Laboratory 
experiments  (Kipp  and  Grady,  1989). 
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FIGURE  5:  Experimental  results  for  SiC  for  plate  impact  velocities  as  indicated  using  the 
experimental  configuration  of  Figure  8  (Kipp  and  Grady,  1989). 


FIGURE  6:  Experimental  results  for  B4C  for  plate  impact  velocities  as  indicated  using  the 
experimental  configuration  of  Figure  8  (Kipp  and  Grady,  1989).  Note  the 
dispersive  unload  behavior  and  absence  of  apparent  strength  on  unload. 
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ti  =  {1  -  V/V0} 


FIGURE  7:  Least  squares  fit  to  Hugoniot  data  for  SiC  (LANL,  1980)  using  equation  5  for  ti 
values  greater  that  0.15  (solid  line).  The  dashed  line  represent  elastic  perfectly- 
plastic  behavior. 


TI  =  {1  •  V/V0} 


FIGURE  8:  Least  squares  fit  to  experimental  Hugoniot  data  for  B4C  (Gust  et  al,  1971)  using 
equation  5  for  T|  values  greater  than  0.15  (solid  line).  The  dashed  lines  represent 
elastic  perfectly-plastic  behavior. 


FIGURE  9: 


Comparison  of  PC-WONDY  model  predictions  and  experimental  data  (Kipp  and 
Grady,  1989)  for  SiC  using  an  elastic  perfectly-plastic  model;  (a)  for  plate  impact 
velocity  of  2100  m/s  and  (b)  for  plate  impact  velocity  of  1542  m/s.  This  model 
shows  clearly  an  absence  of  accurate  yield  function,  the  reflection  of  the  LiF 
window  and  the  presence  of  an  anomalous  wave  at  approximately  1.7  ps. 
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FIGURE  10:  Model  predictions  for  2100  m/s  plate  impact  experiments  with  SiC  with  an 
increase  in  the  strain  rate  hardening  term,  increases  the  apparent  HEL,  increases 
the  magnitude  of  the  elastic  unload  and  causes  the  unload  wave  to  arrive  earlier. 
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FIGURE  12: 


Comparison  of  PC-WONDY  model  predictions  with  experimental  results  for  SiC 
with  plate  impact  velocity  of  2100  m/s  using  the  parameters  C,  =  80,  C,  =  100 

“  “  ~  a  =  0,7x10  kS/J>  and  C3  =  0.008  with  the  unload  isentrope  having  30 

kbar  offset  at  zero  strain.  As  shown,  these  results  compare  favorably  with  the 
experimental  data  (Kipp  and  Grady,  1989). 
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FIGURE  13:  Comparison  of  PC-WONDY  model  predictions  with  experimental  results  for  SiC 
with  plate  impact  velocity  of  1540  m/s  using  the  parameters  Cj  =  80,  C2  =  100, 
n  =  0.1,  a  =  0.7x1c6  kg/J,  and  C3  =  0.008  with  the  unload  isentrope  having  30 
kbar  offset  at  zero  strain.  As  shown,  these  results  compare  favorably  with  the 
experimental  data  (Kipp  and  Grady,  1989)  for  1.2  microseconds  after  plate  impact 
Goading  waves),  but  does  not  model  the  unload  adequately. 
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FIGURE  15:  Calculated  pressure  distribution  versus  position  for  SiC  with  a  plate  impact 
velocity  of  2100  m/s  when  the  pressure  reaches  its  absolute  minimum  in  the  target 
plate. 
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FIGURE  16:  Calculated  pressure  distribution  for  SiC  with  a  plate  impact  velocity  of  2100  m/s 
when  the  lithium  fluoride  window  has  been  removed.  Note  the  signature  of  a 
spall  zone  at  a  position  6  to  7  mm  from  the  flyer/target  plate  interface. 
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FIGURE  17:  Model  predictions  with  an  increase  in  the  damage  coefficient  a  using  the 
parameters;  Cj  =  80,  C2  =  100,  n  =  0.1,  and  Cj  =  0.008  for  SiC  with  plate  impact 
velocity  of  2100  m/s  reduces  the  elastic  unload  and  shows  a  more  dispersive 
plastic  unload  wave. 
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FIGURE  18:  Comparison  of  PC-WONDY  model  predictions  with  experimental  data  for  SiC 
with  plate  impact  velocity  of  1540  m/s  and  with  increased  damage  accumulation. 
The  damage  was  increased  from  6  -  8  percent  (see  figure  15)  to  40  percent.  This 
produces  good  agreement  to  the  experimental  data  shown  with  the  heavy  line 
(Kipp  and  Grady,  1989),  for  times  greater  than  1.3  |is  after  plate  impact. 
Suggestion  is  that  after  the  loading  phase,  damage  continues  to  accumulate, 
reaching  a  level  of  about  40  percent  by  the  time  unloading  begins. 
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r\  =  {1  •  V/VJ 


FIGURE  19:  Least  squares  fit  of  Hugoniot  data  for  B4C  (Gust  et  al,  1971)  using  equation  5  for 
values  of  T|  greater  than  0.15  (solid  line).  The  observed  At|  at  constant  stress  for 
T]  =  .055  to  .065  is  indicative  of  possible  loss  of  porosity  or  phase  change  in  the 
ceramic.  The  shape  of  the  underlying  hydrostat  is  therefore  not  known. 
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FIGURE  20:  Comparison  of  PC-WONDY  model  predictions  with  experimental  data  for  B4C 
using  an  elastic  perfectly-plastic  model  for  plate  impact  velocity  of  1550  m/s. 
Experimental  results  shown  with  the  heavy  line  (Kipp  and  Grady,  1989). 
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(d)  a  =  3.0xl0'5. 


FIGURE  21:  Model  predictions  of  velocity  history  curves  for  B4C  for  plate  impact  velocity  of 
1550  m/s  with  variation  in  the  damage  sensitivity.  Increases  in  a  delays  the 
arrival  time  of  the  plastic  shock.  The  elastic  unloading  wave  begins  to  catch  up 
to  the  plastic  load,  overtaking  it  in  (d). 
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FIGURE  22:  Comparison  between  model  predictions  and  measured  velocity  histories  for  B4C, 
using  a  =  3.4x(10)*5  kg/J,  Ct  =  116,  and  no  work  hardening  or  strain  rate 
hardening  adequately  fits  the  loading  waves  for  plate  impact  velocity  of  1550  m/s. 
Heavier  curve  is  experimental  data  (Kipp  and  Grady,  1989). 
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FIGURE  23:  Comparison  between  model  predictions  and  measured  velocity  histories  for  B4C, 
using  o  =  3.4x(10)'5  kg/J,  =  116,  and  no  work  hardening  or  strain  rate 
hardening  adequately  fits  the  loading  waves  for  plate  impact  velocity  of  2210  m/s. 
Heavier  curve  is  experimental  data  (Kipp  and  Grady,  1989). 
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APPENDIX  1 


RDA  CERAMIC  FAILURE  MODEL 


TENSILE  DAMAGE  MODEL: 

A  material  is  said  to  be  brittle  under  conditions  in  which  its  ability  to  resist  loads 
decreases  with  increasing  deformation  [1].  In  order  to  solve  high  velocity  impact  problems 
with  a  hydro  code  the  total  stress  is  decomposed  into  hydrostatic  and  deviatoric  components. 
To  characterize  the  brittle  nature  of  ceramics  it  is  necessary  to  modify  the  usual  linear  elastic 
relationships  for  stress  and  strain.  The  expressions  for  deviatoric  and  volumetric  stress  are 
therefore  written  in  the  following  form, 


P  =  3K  (1  -  D)  <j> 

0-1) 

Sjj  =  2G  (1  -  D)  ey 

G-2) 

where  P  is  the  pressure,  Sy  and  ey  are  the  deviatoric  components  of  stress  and  strain,  respec 
tively,  <t>  is  the  volumetric  strain,  K  and  G  are  the  undamaged  bulk  and  shear  moduli  of  the 
material,  and  D  is  a  state  variable  describing  the  internal  damage  in  the  material. 

The  underlying  concept  of  this  model  is  that  of  a  degraded  bulk  modulus  K,  and  a 
degraded  Poisson’s  ratio  v,  (Budiansky  and  O’Connell  [6]),  where 


K  =  K  (1  -  16/9  f(v)  Cd) 


f(v)  =(1  -  v2)  /  (1  -2v) 


(1-3) 

0-4) 


v  =  v  (1  -  16/9  Cd)  (1-5) 

Here  v  is  the  undamaged  Poisson’s  ratio  and  Cd  is  the  crack  density,  or  void  volume  of  the 
material. 

Grady  and  Kipp  [3]  proposed  that  the  number  of  flaws  per  unit  volume  at  a  given 
pressure  could  be  described  statistically  using  a  Weibull  distribution 
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N  =  k  (P/3K)m 


(1-6) 


An  expression  for  fragment  diameter  as  a  function  of  strain  rate  was 
also  derived  by  Grady  [4] 


d  = 


(V20  Klc 


V/3 


l  pc*  J 


(1-7) 


Taylor  et  al  [5]  replaced  the  quantity  <J>  with  <}>  the  maximum  strain  at  failure.  They  also 
postulated  that  the  characteristic  flaw  size  was  proportional  to  the  radius  of  the  fragment  size, 
i.e.,  a  ~  l/2d.  Then,  defining  the  average  volume  of  each  flaw  as  Pa3,  where  P  is  an 
unknown  proportionality  constant,  an  expression  for  crack  density,  independent  of  Budiansky 
and  O’Connell  was  given  as 

Cd  =  N  pa3  (1-8) 

or  expanding, 

5  k  rKlc 

C<)  2  (3K)m  =  IjjC 

where, 

P  =  tensile  pressure 
k,m  =  Weibull  parameters 
Klc  =  fracture  toughness 
p  =  density 

C  =  sound  speed 

*mflT  =  maximum  strain  rate  at  failure 
Inspection  of  (1-3)  indicates  that  a  convenient  expression  for  damage  is 

D  =  16/9  f(v)  Cd  (1-10) 
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Arbitrarily  setting  a  limiting  value  of  crack  density,  Cd,  to  9/16,  it  follows  that  the  damage 
varies  smoothly  from  0  to  1  as  the  void  fraction  varies  from  0  to  9/16,  or  to  slightly  more 
than  one-half  its  total  volume. 


In  order  to  describe  the  state  of  the  material,  for  a  given  set  of  volumetric  and 
deviatoric  strain  rates,  the  quantities  Pd.  D,  P,  and  Sy  must  be  solved  for.  These  are  time 
integrated  quantities  and  therefore  must  be  put  into  rate  form.  Taking  the  time  derivatives  of 
(1-9),  (1-10),  (1-1),  and  (1-2),  holding  the  strain  rates  constant,  the  final  form  of  the  rate 
expressions  given  by  Taylor  et  al  is, 


—  Lan-l 


^  0 


-2 


max 


D  = 


16 


[f,  (v)C„  +  C^v)  ] 


g  |  *1  (V)Cd  +  Cdfj 

P  =  3K[(1  -  D)  4>  -  0D] 

Sy  =  2G[(1  -  D)  'ey  -  'ey  D] 


(Ml) 

d-12) 

d-13) 

(1-14) 


We  immediately  make  a  substitution  in  1-11,  replacing  with  0  as  was  in  the  original 
formulation  of  Grady  and  Kipp.  Second,  in  order  for  these  rate  equations  to  be  well  behaved 
as  the  material  cycles  from  tension  to  compression  and  back  an  important  modification  must 
be  made.  Both  the  strain  terms,  0  and  ey,  in  equations  1-13,  and  1-14  must  be  regarded  as 
the  integrated  strains  for  a  given  tensile  cycle,  that  is,  strains  that  resulted  in  the  presence  of 
damage  accumulation.  Once  the  material  unloads  or  goes  into  compression  other  stress  rate 
expressions1  are  used  leading  to  solutions  in  stress  and  strain  space  inconsistent  with  what 
would  result  from  equations  1-13  and  1-14  had  they  been  used.  Referring  to  Figure  1-la  and 
1-lb  we  note  that  residual  strains  can  exist  due  to  plasticity  in  both  tension  and  compression. 
For  example  a  condition  of  zero  stress  can  exist  with  finite  strain.  We  can  see  by  inspection 
of  1-13  and  1-14  that  the  presence  of  this  residual  strain  coupled  with  a  high  damage  rate  can 


lrThe  material  is  assumed  to  unload  in  tension,  and  load  and  unload  in  compression  linear 
elastically.  Plastic  flow  in  compression  is  governed  by  the  von  Mises  incremental  theory  of 
plasticity. 
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I 

produce  a  significant  negative  stress  rate  even  though  the  imposed  strain  rate  is  positive.  It  is 
this  inconsistency  that  leads  to  accumulation  of  negative  energy  i.e.  products  of  negative 
|  stress  rates  with  positive  strain  rates  and  must  be  corrected. 

Fortunately  the  process  of  correcting  this  situation  is  simple  and  straight  forward.  Any 
time  the  material  enters  a  new  tensile  cycle  the  residual  strain  is  ’trapped’.  All  subsequent 
strains  for  this  tensile  cycle  are  corrected  by  this  amount  The  strains  used  in  the  rate 
equations  are 

*1  =  V  <l>r  (1-17) 

where 

<|)  i  =  local  strain 

4>g  =  global  (solution  strain) 

<|>r  =  residual  (’trapped’)  strain 
|  The  stress  rate  equations  now  become 

P  =  3K[(1  -  D)  <j)v  -  4>V*D] 

Sy  =  2G[(1  -  D)  ey  -  ejj*  D] 

where  star  ’*’  denotes  the  corrected,  or  local  strains.  Now,  when  the  material  goes  through  a 
compression  -  tension  transition,  stress  time  histories  are  smooth  and  well  behaved. 

Numerical  stability  is  regained. 


(1-18) 

(1-19) 
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COMPRESSIVE  DAMAGE  MODEL: 


The  compressive  damage  model  is  a  plasticity  based  engineering  model.  Damage  is 
included  as  a  softening  term  in  Johnson  and  Cook’s  expression  for  the  yield  surface  of  a 
ductile  material  [7], 

Y  =  ((Cj  +  q  £pn)  (1+  C3  In  e*)  (1  -  T*m)  +  C4  P)  (1  -  D)  (1-20) 


where  Cj  is  the  static  yield  strength,  C2  and  n  are  the  strain  hardening  parameters,  C3  is  the 
viscoplastic  (strain  rate)  parameter,  T*  and  m  are  the  thermal  softening  terms,  c4  is  the 
confining  pressure  parameter,  and  (1-D)  is  the  damage  induced  softening  term. 


The  damage  evolution  is  related  to  the  plastic  shear  work  rate  in  compression  and 
allows  for  a  single  damage  sensitivity  parameter,  a,  to  characterize  a  given  materials  damage 
dependence  on  plastic  work  rate. 


h  «Wb_ 
"(1  -  D) 


d-21) 


Here  D  is  the  same  scalar  damage  state  variable  as  was  used  in  tension.  Hence  the  integrated 
value  for  damage  during  a  compressive  cycle  becomes  the  initial  value  for  damage  in  a 
subsequent  tensile  cycle  and  vice  versa.  Note  that  any  damage  evolved  in  tension  will  affect 
the  shear  strength  of  the  material  in  compression  (eq.  1-20)  and  that  damage  evolved  in 
compression  will  affect  the  stiffness  in  tension  (eq.  1-1,  1-2).  In  either  case,  when  the 
damage  reaches  one  the  material  loses  all  ability  to  support  either  a  tensile  bulk  or  shear  load, 
or  a  compressive  shear  load.  Program  logic  however,  still  allows  the  material  to  sustain  pure 
hydrostatic  compressive  loads. 
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1-la.  Typical  load  cycle  for  hydrostatic  component  of  stress 
for  ceramic  model.  Damage  evolution  is  reflected  in 
modulus  change.  Unload  path  is  arbitrary  and  is  used 
to  ascribe  dilatancy. 


1-lb.  Typical  load  cycle  for  deviatoric  component  of  stress 
for  ceramic  model.  Damage  accumulation  affects  both 
shear  modulus  in  tension  and  yield  surface  in  compression. 


APPENDIX  2 
FORMULARY 


Bulk  sound  speed: 

C02  =  CL2  -  4/3Cs2 

Poisson’s  ratio: 

CL2  -  2C  2 

V  =  L.  s 

2(Cl2  -  Cs2) 

Isentropic  (cold  compression)  pressure: 

p  =  P0usup 

Assuming 

us  =  C0  +  sup 

then 

(1-STl)2 

Shock  impedance: 

Zl  =  Po^L  (for  cei  ~  ic  with  strength) 


c  _V3(l-v) 
L  (1+v) 


u 


s 


Zl  =  p0us  (for  LiF  with  no  strength) 


(2-1) 

(2-2) 

(2-3) 

(2-4) 

(2-5) 

(2-6) 

(2-7) 

(2-8) 

(2-9) 
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